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Abstract 

We show that any second order linear ordinary diffrential equation 
with constant coefficients (including the damped and undumped har- 
monic oscillator equation) admits an exact discretization, i.e., there 
exists a difference equation whose solutions exactly coincide with solu- 
tions of the corresponding differential equation evaluated at a discrete 
sequence of points (a lattice). Such exact discretization is found for 
an arbitrary lattice spacing. 

1 Introduction 

The motivation for writing this paper is an observation that small and ap- 
parently not very important changes in the discretization of a differential 
equation lead to difference equations with completely different properties. 
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By the discretization we mean a simulation of the differential equation by a 
difference equation p. 

In this paper we consider the damped harmonic oscillator equation 

X + 27X + u^x = . (1) 

where x = x{t) and the dot means the t-derivative. This is a linear equation 
and its general solution is well known. Therefore discretization procedures 
are not so important (but sometimes are applied, see [2]). However, this 
example allows us to show and illustrate some more general ideas. 

The most natural discretization, known as the Euler method fAppendix lB} 
compare [H E]) consists in replacing x by i by the difference ratio 
{xn+i — Xn)ls, X by the difference ratio of difference ratios, i.e., 

1 / ■^n+l -^n+l ■^n\ 2x^+1 ~l~ Xn 

X ^ -[ = , 2 

and so on. This possibility is not unique. We can replace, for instance, x by 
Xn+ij or X by ( )/e, or X by {xn+i - 2x„ + Xn-i)/e^. Actually the 

last formula, due to its symmetry, seems to be more natural than Q (and it 
works better indeed, see Section |21). 

In any case we demand that the continuum limit, i.e., 

Xn = x{tn) , tn = en , e , (3) 

applied to any discretization of a differential equation yields this differential 
equation. The continuum limit consists in replacing Xn by x{tn) = x{t) 
and the neighbouring values are computed from the Taylor expansion of the 
function x(t) at t = t„ 

Xn+k = x{tn + ke) = x{tn) + x{tn)ke + ^x{tn)k'^e^ + ... 

Substituting these expansions into the difference equation and leaving only 
the leading term we should obtain the considered differential equation. 

In this paper we compare various discretizations of the damped (and 
undamped) harmonic oscillator equation. The main result of the present 
paper consists in finding the exact discretization of the damped harmonic 
oscillator equation (Q). By exact discretization we mean that a;„ = x(t„) 
holds for any e and not only in the limit 0. 
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2 Simplest discretizations of the harmonic os- 
cillator 



Let us consider the following three discrete equations 



£2 



£2 



+ Xn-l = , (4) 



-+Xn = 0, (5) 



Xn+l 2Xj2 -|- Xn- 



- + Xn+l = , (6) 



where £ is a constant. The continuum limit © yields, in any of these cases, 
the harmonic oscillator equation 

X + X = . (7) 

To fix our attention, in this paper we consider only solutions correspond- 
ing to the initial conditions a:(0) = 0, i;(0) = 1. The initial data for the 
discretizations are chosen in the simplest form: we assume that Xq and Xi 
belong to the graph of the exact continuous solution. 

For small t„ and small e the discrete solutions of any of these equations 
approximate the corresponding continuous solution quite well (see Fig. 
However, the global behaviours of the solutions (still for small £!) are different 
(see Fig. |21). The solution of the equation (0) vanishes at t — *• cxo while the 
solution of (jD) oscillates with rapidly incrising amplitude (all black points 
are outside the range of Fig. |21). Quahtatively, only the discretization © 
resembles the continuous case. However, for very large t the discrete solution 
becomes increasingly different from the exact continuous solution even in the 
case (0) (compare Fig. |21and Fig. El). 

The natural question arises how to find a discretization which is the best 
as far as global properties of solutions are concerned. 

In this paper we will show how to find the "exact" discretization of the 
damped harmonic oscillator equation. In particular, we will present the 
discretization of (|7j) which is better than (0) and, in fact, seems to be the best 
possible. We begin, however, with a very simple example which illustrates 
the general idea of this paper quite well. 
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Figure 1: Simplest discretizations of the harmonic oscillator equation for small t 
and e = 0.02. Black points: solution of Eq. Q), white points: Eq. Q, grey points: 
Eq. ©, thin line: exact continuous solution. 
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Figure 2: Simplest discretizations of the harmonic oscihator equation for large t 
and e = 0.02. Black points: solution of Eq. Q), white points: Eq. Q, grey points: 
Eq. ©, thin line: exact continuous solution. 
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3 The exact discretization of the exponential 
growth equation 

We consider the discretization of the equation x = x. Its general solution 
reads 

x{t) = x(0)e* . (8) 
The simplest discretization is given by 

= Xn . (9) 

e 

This discrete equation can be solved immediately. Actually this is just the 
geometric sequence Xn+i = (1 + e)xn- Therefore 

x^ = il + erxo. (10) 

To compare with the continuous case we write (1 + e)" in the form 

(1 + £)" = exp(n ln(l + e)) = exp(/s:t„) , (11) 

where t„ := en and k, := ln(l+e). Thus the solution (|iup can be rewritten 
as 

Xn = xoe^'" . (12) 

Therefore we see that for k 7^ 1 the continuous solution (jH}, evaluated at tn, 
i.e. 

x(t„) = x(0)e*" , (13) 

differs from the corresponding discrete solution (fT^ . One can easily see that 
< K, < 1. Only in the limit e ^ we have k —>■ 1. 

Although the qualitative behaviour of the "naive" simulation is in 
good agreement with the smooth solutions (exponential growth in both cases) 
but quantitatively the discrepancy is very large for t 00 because the 
exponents are different. 

The discretization ^ can be easily improved. Indeed, replacing in the 
formula (fTUI) 1 + e by we obtain that (fTUI) coincides with the exact solution 
(Hn)). This "exact discretization" is given by: 

Xn+l Xn A\ 

ee_i = ' (14) 

or, simply, Xn+i = e^x„. Note that ~ 1 + e (for e ^ 0) and this approxi- 
mation yields the equation ©. 
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4 Discretizations of the harmonic oscillator: 
exact solutions 

The general solution of the harmonic oscillator equation (|7j) is well known 

= x(0) cost + i(0) sin t . (15) 

In Section 121 we compared graphically this solution with the simplest discrete 
simulations: (jH), (0), ©• Now we are going to present exact solutions of 
these discrete equations. 

Because the discrete case is usually less known than the continuous one 
we recall shortly that the simplest approach consists in searching solutions 
in the form x„ = A" (this is an analogue of the ansatz x{t) = exp(A)f:) made 
in the continuous case, for more details see Appendix [X]) . As a result we get 
the characteristic equation for A. 

We illustrate this approach on the example of the equation resulting 
from the Euler method. Substituting x„ = A" we get the following charac- 
teristic equation 

A2-2A+(l+e^) = , (16) 
with solutions Ai = 1 + ie, A2 = 1 — ie. The general solution of (jD) reads 

Xn = CiA^ + C2A2 , 

and, expressing ci, C2 by the initial conditions xo,xi, we have 

(1 + ieY - (1 - ier (1 + ie)(l - ie)"' - (1 - ie)(l + ie)" 

Xn = xi — h xo — . 

2ie 2ie 

Denoting 

l + ie = pe*", (17) 



where p = vT+l^ and a = arctane, we obtain after elementary calculations 
Xn = p" (^xq cos{na) + — — — sin(?T,a)^ . (18) 
It is convenient to denote p = e'^^ and 



and then 

Xn = e " [xq COS coin -\ sm totn)- (20) 



One can check that k > and u < 1 for any e > 0. For e — > we have 
K — > 0, ^ 1. Therefore the discrete solution (j^Uj) is characterized by the 
exponential growth of the envelope amplitude and a smaller frequency of 
oscillations than the corresponding continuous solution (|T5|) . 

A similar situation is in the case (jHl), with only one (but very important) 
difference: instead of the growth we have the exponential decay. The formulas 
()19|) and need only one correction to be valid in this case. Namely, k 
has to be changed to — k. 

The third case, (0), is characterized by p = 1, and, therefore, the ampli- 
tude of the oscillations is constant (this case will be discussed below in more 
detail). 

These results are in perfect agreement with the behaviour of the solutions 
of discrete equations illustrated at Fig. [T] and Fig. |21 

Let us consider the following family of discrete equations (parameterized 
by real parameters p, q): 

\-pXn-i + [l-p- q)Xn + qxn+i = , (21) 

The continuum limit Q applied to (^1]) yields the harmonic oscillator ^ for 
any p, q. The family (j2^ contains all three examples of Section |21 and (for 
p = g = 1/4) the equation resulting from the Gauss-Legendre-Runge-Kutta 
method (see Appendix in}: 

Xn+i - 2 (^^^ + = . (22) 

Substituting x„ = A" into (PT|) we get the following characteristic equation: 

(1 + qe')k^ - (2 + (p + g - \)e^)K + (1 + pe^) = (23) 

We formulate the following problem: find a discrete equation in the family 
mP with the global behaviour of solutions as much similar to the continuous 
case as possible. 

At least two conditions seem to be very natural in order to get a "good" 
discretization of the harmonic oscillator: oscillatory character and a constant 
amplitude of the solutions (i.e., p = 1, k = 0). These conditions can be easily 
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Figure 3: Good discretizations of the harmonic oscillator equation for large t and 
e = 0.02. Black points: exact discretization (|l5|) . white points: Eq. ©, grey 
points: Runge-Kutta scheme (|22|) . thin line: exact continuous solution. 
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Figure 4: Good discretizations of the harmonic oscillator equation for small t and 
e = 0.4. Black points: exact discretization (|l5|) . white points: Eq. ©, grey points: 
Runge-Kutta scheme ()2'2() . thin line: exact continuous solution. 
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expressed in terms of roots (Ai, A2) of the quadratic equation ((221) • First, 
the roots should be imaginary (i.e., A < 0), second, their modulus should be 
equal to 1 i.e., Ai = e*", A2 = e~*°. Therefore 1 + pe'^ = 1 + ge^, i.e., q = p- 
In the case q = p the discriminant A of the quadratic equation is given 
by 

A = -Ae^ + e^{l - 4p) . 

There are two possiblities: if p > 1/4, then A < for any e ^ 0, and if 
p < 1/4, then A < for sufficiently small e, namely e"^ < 4(1 — 4p)~^. In any 
case, these requirements are not very restrictive and we obtained p-family of 
good discretizations of the harmonic oscilltor. If Ai = e*" and A2 = e~*", 
then the solution of (PT|) is given by 

xi — xq cos a . / \ /„ ,x 

Xn = xo cositnUJj H sm( t^cj) (241 

sma 

where u = a/e, i.e., 

e - — 



1 

uj = - arctan 

€ 



(25) 



Note that the formula (|24j) is invariant with respect to the transformation 
a — > —a which means that we can choose as Ai any of the two roots of (|23p. 

The equation (0) is a special case of (PT|) for g = p = 0. As we have seen 
in Section |2l for small e this discretization simulates the harmonic oscillator 
((Tj) much better than (jl]) or However, for sufficiently large e (namely, 
6 > 2) the properties of this discretization change dramatically. Its generic 
solution grows exponentially without oscillations. 

Expanding (^3)) in the Maclaurin series with respect to e we get 

1 _u 2 3-40p + 240p2 4 

tu ^ 1 e H € + . . . (26) 

24 640 ^ ' 

Therefore the best approximation of ((7j) from among the family ()2H) is char- 
acterized hy p = 1/12: 

/12 - 5e'^\ 

Xn+l - 2 f _|_ g-2 j^n + Xn~l = ■ (27) 

In this case uj ^ 1 + e'^/ASO + ... is closest to the exact value uj = 1. 
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The standard numerical methods give similar results (in all cases pre- 
sented in Appendix IH] the discretization of the second derivative is the sim- 
plest one, the same as described in our Introduction). The corresponding 
discrete equations do not simulate ((Tj) better than the discretizations pre- 
sented in Section 121 



Damped harmonic oscillator and its dis- 
cretization 



Let us consider the damped harmonic oscillator equation ((H). Its general 
solution can be expressed by the roots Ai,A2 of the characteristic equation 
+ 27A + cUq = and the initial data x(0), x{0): 

\ Ai — A2 / \ A2 — Ai / 

In the weakly damped case {ujq > 7 > 0) we have Ai = —7 + iu and 
A2 = —7 — iuj, where u = \JujI — 7^. Then 

x{t) =x(0) e"^*coswt + tu~^(i;(0) +7x(0)) e'^* sinut . (29) 

To obtain some simple discretization of (0) we should replace the first deriva- 
tive and the second derivative by discrete analogues. The results of Section |21 
suggest that the best way to discretize the second derivative is the symmetric 
one, like in Eq. (0). There are at least 3 possibilities for the discretization 
of the first derivative leading to the following simulations of the damped 
harmonic oscillator equation: 

h 27 h UjQXn = . (30) 

^5 + 27 — +UJ^^Xn = ^. (31) 

2Xn ~l~ Xn—1 , ^n+l -^n , 2 n /on\ 

5 h 27 h uJ^Xn = . (32) 

As one could expect, the best simulation is given by the most symmetric 
equation, i.e., Eq. (jHH), see Fig. El 
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Figure 5: Simplest discretizations of the weakly damped harmonic oscillator equa- 
tion (wo = 1, 7 = 0.1) for small t and e = 0.3. Black points: Eq. (|30j) . white points: 
Eq. grey points: Eq. ((SH), thin line: exact continuous solution. 
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6 The exact discretization of the damped har- 
monic oscillator equation 

In order to find the exact discretization of (Q) we consider tlie general linear 
discrete equation of second order 

Xn+2 = 2AXn+l + Bx^ • (33) 

The general solution of has the following closed form (compare Ap- 
pendix 

_ xo(AiA^-A2A'i')+xi(A?-A^) 
"""" A1-A2 ^ ' 

where Ai, A2 are roots of the characteristic equation A^ — 2AK — B = 0, i.e., 

Ai = A + VA^ + B , A2 = A- VA^ + B . (35) 

The formula (jH^ is valid for Ai 7^ A2, which is equivalent to A"^ + B ^ 0. If 
the eigenvalues coincide (A2 = Ai, B = —A"^) we have Ai = A and 

Xn = il- n)A1xo + nA^^^xi . (36) 

Is it possible to identify Xn given by |^|) with x{tn) where x{t) is given 
by ^? 

Yes! It is sufficient to express in an appropriate way Ai and A2 by Ai and 
A2 and also the initial conditions x{0),x{0) by xo,xi. It is quite surprising 
that the above identification can be done for any e. 

The crucial point consists in setting 

A^ = exp(nlnAfc)) = exp(t„Afc) , (37) 

where, as usual, tn '■= ne. It means that 

Afc:=£-MnAfc, (38) 

(note that for imaginary A^, say A^, = p^e*"*, we have InA^ = Inp^ + iak)- 
Then fj34|) assumes the form 
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Comparing (PS|) with (jSHl) we get x„ = provided that 

x(0) = xo , x(0) = ..A, _ ,sx, • (40) 



The degenerate case, Ai = A2 (which is equivalent to Ai = A2) can be 
considered analogically (compare AppendixEJ. The formula flHUj) is obtained 
from ()34p in the limit A2 Ai. Therefore all formulas for the degenerate 
case can be derived simply by taking the limit A2 — ^ Ai. 

Thus we have a one-to-one correspondence between second order differen- 
tial equations with constant coefficients and second order discrete equations 
with constant coefficients. This correspondence, referred to as the exact dis- 
cretization, is induced by the relation (jHHj) between the eigenvalues of the 
associated characteristic equations. 

The damped harmonic oscillator corresponds to the discrete equation 
fjHHjl such that 



2A = e-'^ (e"VT'-^o + e'^V^^A , B = -e'^"^ . (41) 



In the case of the weakly damped harmonic oscillator {ujq > 7 > 0) the 
exact discretization is given by 

A = e-'"' cos(ew) , B = -g-^^"^ , (42) 



where a; := yt^Q — 7^. In other words, the exact discretization of ((T)) reads 

Xn+2 - 2e~^'^ cos(cje)x„+i + e"^^^x„ = . (43) 
The initial data are related as follows (see pUj) ): 

. xiwe^^ - (7 sm{uje) + lo co?,{uje))xo 
x{0) = xo , a:(0) = , 

(44) 

/ , , sin(cij£:) / sm(uje) , \ , \ 

= ix{0) — ^ + (7 — ^ + cos(tue) j x{0) j e"^^ . 

Fig. El and Fig. [7| compare our exact discretization with two other good dis- 
cretizations of the weakly damped harmonic oscillator. The exact discretiza- 
tion is really exact, i.e., the discrete points belong to the graph of the exact 
continuous solution (for any e and any n). Similarly as in the undumped 
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Figure 6: Good discretizations of the weakly damped harmonic oscillator equation 
(wo = 1) 7 = 0.1) for small t and e = 0.2. Black points: exact discretization ()45() . 
white points: Eq. grey points: Runge-Kutta scheme H22|l . thin line: exact 

continuous solution. 
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Figure 7: Good discretizations of the weakly damped harmonic oscillator equation 
(wo = 1, 7 = 0.1) for large t and e = 0.2. Black points: exact discretization ()45() . 
white points: Eq. (|31|) . grey points: Runge-Kutta scheme H22|l . thin line: exact 
continuous solution. 
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case, the fully symmetric simple discretization is better than the equa- 
tion resulting from GLRK method. 

The exact discretization of the harmonic oscillator equation x + x = is 
a special case of (jlHj) and is given by 

Xn+2 - 2(cos£)x„+i + x„ = . (45) 

It is easy to verify that the formula (|45|) can be rewritten as 

2x^+1 ~l~ X^ 



(2sin(£/2))2 



+ Xn+l = , (46) 



which reminds the "symmetric" version of Euler's discretization scheme (see 
(121) and (0)) but e appearing in the discretization of the second derivative is 
replaced by 2sin(£:/2). For small e we have 2sin(£:/2) ^ e. 

The comparison of the exact discretization with three other discrete 
equations simulating the harmonic oscillator is done in Fig. |3J and in Fig. |3] 
We point out that the considered simulations are very good indeed (although 
t in Fig 121 is very large) but they cannot be better than the exact discretiza- 
tion. The discretization pTj) is also excellent. The coefficient by —2xn in 
Eq. ^ 

12 -5^2 1 1 

T2T7^"^-r 

approximates cose up to 4th order. Actually, for the choice of parameteres 
made in Fig. El and Fig. HI the discretization pTjl practically cannot be dis- 
cerned from the exact one. 



7 Conclusions 

In this paper we have shown that for linear ordinary differential equations 
of second order with constant coefficients there exists a discretization which 
is "exact" and simulates properly all features of the differential equation. 
The solutions of this discrete equation exactly coincide with solutions of the 
corresponding differential equation evaluated at a discrete lattice. Such exact 
discretization can be found for an arbitrary lattice spacing e. 

Therefore we conclude that in this case differential and difference equa- 
tions are in one-to-one correspondence: to any linear differential equation 
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with constant coefficients there corresponds a difference equation which we 
call the exact discretization. 

Analogical considerations can be made for linear ordinary differential 
equations (with constant coefficients) of any order (the details will be pre- 
sented elsewhere). 

We point out that to achieve our goal we had to assume an essential de- 
pendence of the discretization on the considered equation, in contrast to the 
standard numerical approach to ordinary differential equations where prac- 
tically no assumptions are imposed on the considered system (i.e., universal 
methods, applicable for any equation, are considered) 

In the last years one can observe the development of the numerical meth- 
ods which are applicable for some clases of equations (e.g., admitting Hamil- 
tonian formulation) but are much more powerful (especially when global or 
qualitative aspects are considered) |31 IH] . 

"Recent years have seen a shift in paradigm away from classical consider- 
ations which motivated the construction of numerical methods for ordinary 
differential equations. Traditionally the focus has been on the stability of dif- 
ference schemes for dissipative systems on compact time intervals. Modern 
research is instead shifting in emphasis towards the preservation of invariants 
and the reproduction of correct qualitative features" 171 . 

In our paper we have a kind of extremal situation: the method is ap- 
plicable for a very narrow class of equations but as a result we obtain the 
discretization which seems to be surprisingly good. 

Similar situation occurs for the integrable (soliton) nonlinear systems. It 
is believed (and proved for a very large class of equations) that integrable 
equations admit integrable discretizations which preserve the unique features 
of these equations (infinite number of conservation laws, solitons, transfor- 
mations generating exphcit solutions etc.) [HIE]- 

The exact discretization considered in this paper is the best possible simu- 
lation of a differential equation. Linear ordinary differential equations always 
admit the unique exact discretization. An open problem is to find such dis- 
cretization for some other classes of differential equations. 
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Appendix 



A Linear difference equations with constant 
coefficients 

We recall a method of solving difference equations with constant coeficients. 
It consists in representing the equation in the form of a matrix equation of 
the first order. The general linear discrete equation of the second order 

Xn+2 = 2Ax„+i + BXn , (47) 

can be rewritten in the matrix form as follows 

Vn+i = Myn , (48) 

where 



The general solution of (j48|) has, obviously, the following form: 

Vn = M^y, (50) 

and the solution of a difference equation is reduced to the purely algebraic 
problem of computing powers of a given matrix. 

The same procedure can be applied for any linear difference equation 
with constant coefficients. If the difference equation is of m-th order, then 
to obtain the equation (P5|) we define 

Un ■ (•^n+m,) •^m+m—li • • • i -^n+li -^n) (^-^) 

where the superscript ^ means the transposition. 

The power M" can be easily computed in the generic case in which the 
matrix M can be diagonalized, i.e., represented in the form 



-1 



M = NDN 

where Z) is a diagonal matrix. Then, obviously, M" = ND^N~^. The 
diagonalization is possible whenever the matrix M has exactly m linearly 
independent eigenvectors (in particular, if the characteristic equation 
has m pairwise different roots). Then the columns of the matrix are 
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just the eigenvectors of M , and the diagonal diagonal coefficients of D are 
eigenvalues of M. 

The characteristic equation (det(M — XI) = 0) for m = 2 (i.e., for (jUj)) 
has the form 

- 2A\ -B = . (52) 

Its roots will be denoted by Ai, A2 (see (jHSl))- If ^1 7^ ^2, then the diagonal- 
ization procedure yields 

where the columns of are the eigenvectors of M, i.e., 

^-(\'\')- (54) 
Therefore 

and performing the multiplication we get ()34p. 

The case of multiple eigenvalues of M is technically more complicated. 
In order to compute M" we can, for instance, transform M to the Jordan 
canonical form (see, for instance, |10j). Here we suggest a method which is 
very efficient for 2x2 matrices. By the Cayley-Hamilton theorem ( JU]) any 
matrix satisfies its characteristic equation. In the case of ()47|) it means that 
= 2AXM + B. In the case of the double root {B = —A^) one can easily 
prove by induction 

M" = (1 - n)A"M + nA"-^ . (56) 

Substituting it to (pH|) we get immediately (|3ti|l. 



B Numerical methods for ordinary differen- 
tial equations 

In this short note we give basic informations about some numerical methods 
for ordinary differential equations and we appply them to case of harmonic 
oscillator equation ((Tj). 
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A system of linear ordinary differential equations (of any order) can al- 
ways be rewritten as a single matrix equation of the first order: 

y = Sy, (57) 

where the unknown ?/ is a vector and S" is a given matrix (in general t- 
dependent). Numerical methods are almost always (see j^) constructed for 
a large class of ordinary differential equations (including nonlinear ones): 

y = f{t,y). (58) 

We denote by ?/„ a numerical approximant to the exact solution ?/(t„). 

Euler's method 

yk+i=yk + ef{tk,yk) (59) 
In this case the discretization of x + a: = is given by (jlj). 

Modified Euler's methods 

Vk+i = yk + ef{tk + ^e, yk + ^^/(tfc, yk)) (60) 

yk+i = yk + ^^ifitk, yk) + f{tk + i/fc + £f{tk, yk))) (61) 
Both methods lead to the following discretization of a; + a; = 0: 

5 VXn + -e Xn-l = . (62) 

The roots of the characteristic equation are imaginary and 



|Ai| = |A2| = Wl + -. (63) 
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1-stage Gauss-Legendre-Runge-Kutta method 

Vk+i = yk + efitk + 2 ^ ^ ^ 

The application of this numerical integration scheme yields the following 
discretization of the damped harmonic oscillator equation: 

^ ^ 2i ^ 4 " • ^ ^ 

In the case 7 = 0, a;o = 1 (i.e., x + x — 0), we have 

2 + ie . 2 — i£ , , 

A. = ^. A, = — . (66) 

Adams-Bashforth extrapolation formula 

k 

Vn+l = 2/n + £ X) bkjf{tn-j, Vn-j) (67) 

where hkj are specially chosen real numbers. In particular: 610 = 3/2, fen = 
-1/2, 620 = 23/12, 621 = -4/3, 622 = 5/12. 
In the case A; = 1 we obtain (for x -\-x — Qi): 

-2 + 1 -0 (68) 

and the characteristic equation reads 

- 2A^ + (1 + -s' ] - -e'A + -e' 



This is an equation of the 4th order (with no real roots for s ^ 0). 
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